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ABSTRACT 


The  problem  of  numerical  integration  of  the  parabolic  differential 
equations  of  viscous  shock-layer  theory  are  studied.  Three  methods  are 
considered;  the  standard  explicit  and  implicit  schemes,  as  well  as  the 
explicit  scheme  of  DuFort  and  Frankel.  The  relative  merit  and  efficiency 
of  these  methods  are  discussed  in  terms  of  numerical  solutions  carried  out 
for  the  shock  layer  on  a  cone. 

Because  of  the  singularity  in  the  governing  equation,  the  mesh  size 
associated  with  the  difference  approximation  must  be  progressively  refined 
as  the  leading-edge  is  approached.  This,  in  practice,  places  a  severe  handi¬ 
cap  on  all  difference  methods,  particularly  the  explicit  scheme  which  is 
already  restricted  by  the  requirement  of  stability.  The  examples  worked 
out  here  for  flows  over  cones  confirm  that  the  implicit  and  DuFort-Frankel 
explicit  schemes,  which  are  not  subject  to  the  stability  requirement,  are 
more  effective.  Although  the  difference  between  the  two  latter  schemes  is 
not  large  in  this  example,  the  discussion  reveals  that  the  implicit  scheme, 
which  is  subject  to  neither  the  stability  nor  the  convergence  requirement, 
is  much  more  efficient  and  accurate  in  problems  involving  non-vanishing 
pressure  gradients. 
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I.  INTRODUCTION 


There  are  many  methods  for  obtaining  numerical  solutions  to  nonlinear 
partial  differential  equations  of  hyperbolic  or  parabolic  type.  The  following 
study  concerns  the  application  of  numerical  methods  to  the  study  of  the  low- 
density  hypersonic  flow  based  on  the  mathematical  formulation  of  Ref.  1,  in 
which  the  system  of  governing  equations  is  reduced  to  one  of  the  parabolic 
type,  similar  to  those  of  the  boundary-layer  theory. 

In  preparation  for  a  more  general  attack,  the  present  investigation  is 
made  with  the  objective  of  examining  the  applicability  of  certain  finite-dif¬ 
ference  techniques  to  the  hypersonic  shock-layer  theory.  For  this  purpose, 
analyses  and  discussions  will  be  made  only  for  the  simplest  example,  that 
of  hypersonic  flow  over  a  cone  at  low  Reynolds  number.  In  this  case,  under 
the  assumption  of  a  unit  Prandtl  number,  the  governing  equations  are  reducible 
to  a  single  parabolic  equation,  but  the  solution  is  no  long'er  self-similar  as 
in  the  boundary-layer  theory  because  of  the  modified  outer  boundary  condition. 
Near  the  cone  apex,  an  analytic  development  in  ascending  powers  of  the  dis¬ 
tance  from  the  apex  is  valid;  while  far  downstream,  an  asymptotic,  self¬ 
similar  solution  corresponding  to  the  classical  boundary-layer  theory  exists. 
The  finite-difference  methods  provide  a  valid  transition  between  the  two 
solutions.  The  agreement  with  these  two  solutions  at  small  and  large  dis¬ 
tances  provides  a  check  of  the  accuracy  of  the  finite -difference  methods. 

In  the  method  of  finite  differences,  the  derivatives  are  replaced  by  the 
quotients  of  the  differences,  and  the  partial  differential  equations  are  satisfied 
approximately  at  a  finite  number  of  grid  or  lattice  points.  For  an  initial-value 
problem  pertaining  to  the  parabolic  and  hyperbolic  equations,  solution  by 
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forward  integration  is  possible,  and  the  solution  to  the  difference  equation  can 
be  determined  in  terms  of  previously  obtained  values  at  certain  preceding  sta¬ 
tions  or  grid  line.  By  employing  sufficiently  small  spacing  between  grid  points, 
and  with  the  aid  of  a  high-speed  digital  computer,  one  may  obtain  a  numerical 
solution  with  high  accuracy. 

The  manner  in  which  the  difference  quotients  are  formed  and  the  equations 

(2  31 

approximated  gives  rise  to  a  variety  of  difference  methods.  '  ’  '  The  simplest 

and  most  commonly  known  of  these  is  the  standard  explicit  scheme  by  which 

each  unknown  quantity  is  determined  explicitly  in  terms  of  values  obtained  pre- 

viously  at  a  preceding  station.  '  ’  '  This  method  is,  however,  restricted  by 

the  stability  requirements.  Another  scheme  introduced  by  DuFort  and  Frankel^^’ 

is  also  explicit,  in  that  each  unknown  value  of  the  solution  is  determined  explicitly 

with  the  known  values  at  two  preceding  stations  in  a  specific  manner.  Although 

this  method  is  generally  stable,  the  question  of  convergence  may  become  a  crit- 

(2  31 

ical  problem.  In  an  implicit  scheme,  '  ’  '  each  unknown  value  of  the  solution  is 
not  determined  explicitly  in  terms  of  the  known  quantities  at  the  preceding  grid 
line  but  has  to  be  determined  simultaneously  with  other  unknowns  along  the  same 
grid  line.  The  relative  advantage  in  this  method  is  that  there  is  no  stability  and 
convergence  requirement.  For  brevity,  the  explicit,  DuFort-Frankel  explicit 
and  implicit  schemes  will  be  referred  to  hereafter  as  the  E,  D-F  and  I  schemes, 
respectively.  These  three  schemes,  which  are  most  widely  used,  are  to  be  dis¬ 
cussed.  The  variants  of  these  schemes  as  well  as  other  methods  will  not  be 
pursued  here. 

The  finite -difference  methods  have  been  applied  to  viscous  flow  problems 

lS-9) 

within  the  context  of  the  boundary-layer  theory  by  a  number  of  investigators.  '  ' 
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Wu  used  the  E  scheme,  and  Der  and  Raetr^^^  employed  the  D-F  scheme. 
Kramer  and  Lieberstein^^^  and  Flugge-Lotz  and  co-workers^®’  have  applied 
and  studied  E  and  I  schemes  as  well  as  their  variants  to  a  number  of  boundary- 
layer  problems.  None  has  studied  and  compared  all  three  schemes.  In 
addition,  the  singular  behavior  of  the  differential  and  difference  equations  of 
the  present  study  differs  from  the  singularities  in  the  cited  boundary-layer 
work  because  of  the  different  formulations  used.  It  is  essential  to  observe  in 
this  connection  that,  in  many  boundary -layer  works,  the  momentum -integral 
and  other  equivalent  methods  may  be  sufficient.  But  for  the  present  pro¬ 
blem,  as  well  as  other  hypersonic  rarefied- gas  flow  problems,  the  integral 
methods  (which  also  involve  the  stepwise  numerical  calculations)  are  shown  to 
be  inadequate  because  the  velocity  profiles  do  not  permit  simple  description 
over  a  wide  range  of  Reynolds  number. 

In  the  subsequent  section,  certain  preliminary  remarks  on  the  errors, 
stability  and  convergence  of  the  difference  approximations  will  be  given.  The 
model  of  the  flow  problem  and  its  governing  equations  are  described  in  Section 
3.  In  Section  4,  essential  details  of  the  three  finite  difference  methods  are 
discussed.  In  Section  5,  the  results  obtained  with  these  methods  are  examined 
and  the  relative  merit  of  the  three  schemes  are  evaluated. 
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II.  STABILITY  AND  CONVERGENCE 
In  general,  there  are  two  types  of  errors  which  cause  departures  from 
the  exact  solution  to  the  partial  differential  equation;  truncation  errors  and 
round-off  errors.  The  first  arises  from  replacing  the  derivatives  by  the  dif¬ 
ference  quotients  and  may  be  considered  as  an  error  resulting  from  the  trun¬ 
cation  of  the  Taylor  series.  The  second  results  from  the  rounding-off  in  the 
numerical  computation.  Generally,  the  effect  of  these  errors  may  be  controlled 
by  reducing  the  spacing  between  grids  and  including  more  significant  figures  in 
the  computation. 

In  practice,  however,  the  departure  from  the  exact  solution  due  to  each 

of  these  errors  may  grow  or  decay  as  computation  progresses  in  a  manner  which 

2 

is  a  property  of  the  difference  scheme.  Under  certain  circumstances,  the  rate 

of  growth  of  the  departure  may  become  too  large  to  be  manageable,  that  is,  it 

becomes  unstable.  On  the  other  hand,  there  are  schemes  which  are  inherently 

stable  in  the  sense  that  in  the  course  of  computation  the  departure  due  to  errors 

tends  to  decay.  Many  schemes  are  only  conditionally  stable;  the  E  scheme  is 

2  3 

one  of  the  commonly  known  examples.  '  The  D-F  scheme  is  generally  stable 

but  subject  to  the  problem  of  convergence.  That  is,  as  the  spacing  of  the  mesh 

is  reduced,  the  solution  of  the  difference  equation  can  not  always  be  reduced  to 

2  3 

that  of  the  original  differential  equation.  ’  The  application  of  either  of  the 
E  and  D-F  schemes  is  rather  straightforward,  but  is  handicapped  by  either  the 
stability  or  convergence  requirement  which  imposes  restrictions  on  the  choices 


Error  may  also  result  from  the  inaccurate  description  of  the  initial  and 
boundary  conditions. 
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of  the  incremente  in  the  two  independent  variables.  In  contrast,  whereas  the 
I  scheme  is  generally  free  from  both  of  the  stability  and  convergence  require¬ 
ments.  ^ 

Apparently,  in  many  cases,  these  restrictions  can  be  satisfied  by  suitable 
adjustments  to  the  spacing  used  for  thu  mesh.  However,  for  a  differential 
equation  such  as  the  one  to  be  considered,  because  of  the  singularity  in  the 
difference  equation,  the  stability  and  convergence  requirements  may  become 
stringent  enough  to  make  further  refinements  in  the  mesh  size  economically 
unfeasible.  This  singularity  will  occur  at  the  apex  of  the  cone,  the  leading 
edge  of  a  flat  plate  or  the  stagnation  region  of  a  blunt  body.  In  order  to  obtain 
the  solution  in  the  vicinity  of  the  singularity,  an  expansion  in  ascending  powers 
of  the  distance  from  the  singular  point  has  to  be  developed  so  that  numerical 
integration  may  begin  at  some  downstream  station.  In  practice,  this  pro¬ 
cedure  introduces  an  error  in  the  initial  data  because  of  the  use  of  truncated 
series  in  the  expansion.  The  error  could,  of  course,  be  reduced  by  initiating 
the  numerical  integration  at  a  station  closer  to  the  singularity,  where  higher- 
order  terms  in  the  expansion  are  not  required.  But,  because  of  the  singular 
behavior  of  the  equations,  moving  the  station  upstream  would  reduce  the 
scheme's  accuracy  as  well  as  stability  (if  the  scheme  is  only  conditionally 
stable).  Thus  the  presence  of  the  singularity  leads  inevitably  to  the  use  of  a 
refined  mesh  size  in  the  vicinity  of  the  singularity  and  to  excessively  long 
computation  time. 
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in.  THE  FLOW  MODEL  AND  ITS  GOVERNING  EQUATION 


The  flow  model  considered  pertains  to  a  continuvim  description  of  hyper¬ 
sonic,  rarefied  gas  flows  over  cones.  For  simplicity,  a  perfect  gas  with  constant 
specific  heats  is  assumed.  Without  going  through  the  details  of  the  formulation 
based  on  the  thin  shock-layer  theory  of  Cheng,  ^  the  differential  equation  govern¬ 
ing  the  velocity  under  the  assumption  of  a  linear  viscosity-temperature  relation 
can  be  written  for  a  cone  or  wedge  in  terms  of  a  pair  of  dimensionless  variables 
and  as 


2 


(1) 


The  inner  boundary  condition,  at  =0,  is 


W(0)  -0  (2) 

and  the  outer  boundary  condition  immediately  behind  the  shock,  at  y  =  1,  is 

(3) 

wh.t.  0^  and  K  are  the  free- stream  and  local  velocities 

respectively.  The  index  3/  is  zero  for  a  wedge  and  unity  for  a  cone.  For  the 
present  study,  if  is  therefore  always  taken  as  unity.  The  effects  of  the  gas 
rarefaction  as  well  as  surface  inclination  are  absorbed  in  the  variable 
7  S  which  is  essentially  a  function  of 

local  Reynolds  number.  The  variable  y  is  related  to  the  stream  function  y 
through  y  S.  (n'Zj.  One  advantage  of  using  y  as  an  indepen¬ 

dent  variable  is  that  it  eliminates  the  need  for  determining  the  location  of  the 
outer  edge  of  the  shock  layer,  as  is  evident  from  the  boundary  condition  Eq.  (3). 

After  the  velocity  is  obtained,  the  total  enthalpy  can  in  turn  be  found. 

For  a  unit  Prandtl  number  the  Crocco  relation  gives 
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u 


where  subscripts  oo  and  U/  denote  the  total  enthalpies  pertaining  to  free- 
stream  and  the  surface  conditions,  respectively. 

As  is  characteristic  of  the  equations  governing  boundary  layers  and 
shock  layers  the  coefficient  of  vanishes  at  =  0  which  in  this  case 

is  the  cone  apex.  This  will  give  rise  to  a  singularity  in  any  method  employing 
forward  integration  in  ^  .  In  addition,  because  of  the  presence  of 
the  equation  and  solution  are  also  singular  at  the  surface,  i.  e.  ^  =  0.  The 
latter  singularity  can  be  avoided  by  using  the  combination  of  dependent  and 
independent  variables  U  and  ,  "/W)  >  place  of  IV  and  ( ^  )  )  . 

Subsequent  discussion  will  show  however  that  the  singularity  associated  with 
TV  =  0  turns  out  to  be  unimportant.  In  fact,  corrections  are  not  necessary 
for  the  present  analysis. 

As  pointed  out  in  the  preceding  discussion,  a  development  in  powers  of 
^  will  be  used  to  provide  a  solution  in  the  vicinity  of  the  apex  and  to  furnish 
initial  data  for  numerical  integration  beginning  at  some  downstream  station. 
The  expansion  admissible  by  the  differential  equations  and  boundary  conditions, 
Eqs,  (2)  and  (3),  is  of  the  form 

if/  .  (5) 

or,  alternatively 

“  ■  =  “•  ^  . 

where  the  coefficients  u,,  U,  ,  etc.  are  determined  after  collecting  equal 

powers  of  ,  as 
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IV.  THE  THREE  FINITE -DIFFERENCE  METHODS 
In  the  following,  the  resultant  difference  equation,  in  which  the  difference 
quotients  approximate  the  partial  derivatives,  will  be  described  for  each  of  the 
three  schemes  mentioned.  The  manner  in  which  the  increments  A  7  and  AV 
will  affect  and  control  the  errors,  together  with  the  stability  and  convergence 
characteristics  of  the  schemes,  will  now  be  studied  more  closely.  No  attempt 
has  been  made  to  establish  the  proofs  of  stability  and  convergence  for  these 
schemes,  which  is  clearly  not  the  objective  of  this  report;  however,  the  stab¬ 
ility  and  convergence  for  each  of  the  schemes  are  demonstrated  by  the  solutions 
obtained. 

1,  The  Standard  Explicit  Scheme  (E-Scheme) 

(A)  Difference  Equation  -  For  clarity,  the  reader  is  referred  to  the 
sketch  of  Fig.  la.  In  this  scheme,  the  data  along  a  grid  line,  say  97^  column, 
are  assumed  known;  and  the  partial  derivatives  at  a  point  in  this  column 

*7^Ai?  I  V'-T^Ay )  S'l'e  evaluated  by  the  difference  quotients  which  appear  as 
the  first  terms  in  the  following 

I  / .  . 

^  •  •  •  s  • 


»■»/  Wm.it-/ _ i 

Sy> 


(7) 


(  _  ^rn.nH  ~  ZWm.n  Wm,n~l  I  /a 
T — — ^ 


'W*  (A  f  )* 

The  second  term  on  the  right  of  each  of  the  above  equations  represents  the 
estimate  of  the  truncation  error  based  on  Taylor's  Theorem. 

Upon  using  these  difference  quotients,  Eq.  (1)  is  approximately  satisfied 
at  the  point  (  m  >  Tf  )  and  the  value  of  W  at  the  downstream  neighboring  point 
can  be  evaluated  explicitly  in  terms  of  values  at  the  three  points 
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(mtn+l)  .  ,  7l)  and  (tij  ,  n  -  i) 


w, 


■I. »-/)  7‘/a?V”  »-j)  ( 8 ) 


»i*i,n~  '  }ifV'"i,nH  ”m,n‘ij  '  ^ (A^ 

The  boundary  conditions  in  the  finite-difference  form  will  be  discussed 
separately  later,  as  they  are  applicable  in  the  same  manner  to  all  three  methods. 

(B)  Singularities  and  Trxmcation  Error  -  Assuming  that  the  solution  W 
is  regular  in  both  7  and  ^  ,  application  of  Eq.  (7)  introduces  a  truncation 
error  in  Eq.  (8)  belonging  to  the  order  of  )  J  When  the  step- 

sizes  AS  and  A'V  are  made  successively  smaller,  one  may  anticipate  that  the 
truncation  error  will  diminish  and  the  solution  to  the  difference  equation  Eq. 

(8)  will  converge  to  that  of  the  differential  equation  Eq.  (1).  It  is  important  to 
examine  more  carefully  the  nature  of  the  truncation  error  in  the  neighborhood 
of  ^  =  0,  as  well  as  that  of  ^  =  0  where  the  solution  is  singular. 

The  behavior  of  W  near  the  cone  surface,  i.e.  =  0,  can  be  inferred 


from  the  differential  equation  Eq.  (1)  as 

w  =  o(f 


(9) 

_  Va 


where  CL  and  b  are  functions  of  ^  .  Because  of  the  singular  term  b'V'  , 
the  difference-quotient  approximations  of  and  will  contain  errors 

larger  than  the  order  of  (A1/^)  which,  however,  can  be  corrected.  The 
largest  errors  of  this  kind  occur  at  ^  =  AY*  •  that  is  the  first  grid  point  from 
the  boundary  f  =  0,  Using  Eq.  (9),  one  can  express  the  difference  quotients 
in  terms  of  A  and  b  .  Thus,  at  Y’  s:  AT 
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On  the  other  hand,  Eq.  (9)  also  yields,  at 

Wf^  a- 

It  follows  that,  at  the  worst  situation  {  ) 

ivv  -  JiV  =  o(^^f 

-Wff  =  h(4rz-z-  o(4^f 

The  factors  (  ^  )  and  (  4Jz-’Z  -  ^  )  have  the  values  (+0.  3284)  and 

(— .  0932),  respectively.  The  magnitude  of  is  generally  less  than  unity. 

In  fact,  the  development  given  by  Eq.  (5)  indicates  that  ^  is  considerably 
less  than  1/30.  Therefore,  having  in  mind  an  increment  A'ip'  of  the  order 
of  l/IO,  as  in  most  of  the  present  calculations,  the  truncation  error  resulting 
from  the  singularity  at  ^  =  0  is  numerically  comparable  to,  or  smaller  than, 

.  Consequently,  correction  of  such  an  error  is  not  necessary.  There¬ 
fore,  the  function  will  be  treated  as  being  regular  with  respect  to  all  ^p' 
for  all  three  schemes. 

The  truncation  error  is  affected  more  critically  by  the  singularity  at 
^  =  0  than  that  at  Ap  =  0.  In  fact,  it  is  so  severely  affected  that  the  region 
near  ^  =  0  has  to  be  excluded  from  the  numerical  integration.  To  elucidate 
how  the  truncation  error  depends  on  x  »  may  apply  Eq.  (7)  to  the  dif¬ 
ferential  equation,  Eq.  (1),  retaining  the  remainders  of  the  difference  quotients, 

where  the  functions  ,  and  ^^®  left-hand  side  are  the  dif- 

“X  T 

ference  quotients  given  in  Eq.  (7).  The  right-hand  side  of  Eq.  (11)  represents 
the  estimate  of  the  truncation  errors  resulting  from  the  difference  approximation. 


11 


AF-1285-A-9 


While  %  is  seen  to  appear  in  the  last  two  terms  as  and 
(11),  the  manner  in  which  these  error  terms  will  vary  with  "Z  depend  on 
the  singular  behavior  of  IV  with  respect  to  "x  Applying  the  development 
of  Eq.  (5)  for  small  x  ,  with  the  coefficients  given  by  Eq.  (6),  the  terms 
on  the  left-hand  side  of  Eq.  (11)  are  found  to  be  the  order  of  ^  ,  whereas 

the  three  terms  on  the  right-hand  side  in  succession  are,  respectively,  in 
the  order  of  'ip‘ ’  "Z  ’  •  and  .  The  six  terms 

in  Eq.  (11)  are  therefore  in  the  proportion 


At  a  given  ^  ,  the  last  two  terms  of  Eq.  (12)  become  largest  at  the  smallest 
value  of  ;  at  =  ,*  they  become  of  the  order  and  1, 

respectively.  On  account  of  the  last  of  the  remainders,  the  truncation  error 
of  the  difference  equation  cannot  be  made  arbitrarily  small  as  is  made 

to  vanish.  Fortunately,  the  terms  on  the  right-hand  side  of  Eq.  (11)  associated 
with  are  always  small,  numerically  speaking.  This  fact  can  readily  be 

checked  by  making  use  of  the  coefficients  of  Eq.  (6).  To  represent  more 
appropriately  the  relative  magnitudes  in  a  numerical  sense,  the  last  two  terms 
in  Eq.  (12)  have,  in  fact,  to  be  multiplied  by  a  factor  of  1/30.  Therefore,  from 
a  practical  viewpoint,  the  finite-difference  approximation  may  still  be  applicable 
for  small  % 

(C)  Stability  -  As  pointed  out  before,  the  difficulty  of  the  E  scheme  is 
not  with  the  truncation  error  alone,  rather,  it  is  a  problem  of  stability.  For  a 

♦  ^ 

The  difference  equation  is  not  applied  at  =  0. 


in  Eq. 
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simple  parabolic  equation,  e.  g. 


=  O' 


2^ 


(13) 


with  a  constant  coefficient  (J'  ,  the  criterion  for  the  stability  of  the  E 

scheme  is^ 


Z(r 


AX 


«  L 


The  presence  of  a  lower-order  derivative,  such  as 


(14) 


,  and  to  will 


affect  the  stability  characteristics  of  the  difference  equation  only  through  terms 

of  orders  higher  than  AX  and  ( ;  consequently  they  are  not  important 

^  2  3 

from  a  practical  viewpoint.  Similarly,  the  presence  of  variable  coefficients,  ’ 

as  well  as  the  nonlinearity,  may  also  be  looked  upon  as  higher-order  effects  in 
the  stability  consideration,  provided  that  the  coefficient  O'  be  interpreted 
locally,  and  that  increments  AX  and  be  taken  sufficiently  small. 

It  should  be  noted  at  this  point  that  Eq.  (14)  can  be  shown  to  be  the  nec¬ 
essary  and  sufficient  condition  for  the  stability  of  the  E- scheme  with  a  constant 
coefficient  CT  (see  Refs.  11  and  12).  It  has  also  been  shown  that  Eq.  (14) 

is  a  sufficient  condition  for  stability  of  the  E- scheme  for  a  quasi-linear 

*•  13 

equation 


In  this  connection,  a  criterion  based  on  an  electric  network  analogy  may  be 

of  interest.  By  considering  the  stability  for  the  current  or  voltage  in  an 

14 

electric  network,  Karplus  studies  the  difference  equation  in  the  form 
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> 


(15) 


(Kuo ^  ^ 

. 

^  ^  (K,K'*rK,n)’^  ^(K.>ur  K,J 

^ -  j 

where  a,  b,  c,  d,  and  e  are  the  coefficients  of  the  difference  terms.  It 
has  been  shown  in  Ref.  14  that  if  all  the  coefficients  are  positive,  the 
equation  is  stable;  and  if  some  of  the  coefficients  are  negative,  a  sufficient 
condition  for  stability  is  that  the  algebraic  sum  of  all  the  coefficients  be 
negative.  For  the  simple  parabolic  equation,  Eq.  (13),  application  of  the 
Karplus  condition  to  the  E-scheme  yields  the  criterion  of  Eq.  (14).  The 
Karplus  condition  may  be  quite  useful  in  suggesting  stability  characteristics 
of  new  schemes. 

For  the  stability  of  Eq.  (8)  of  the  E  scheme  the  criterion  Eq.  (14) 
gives  the  condition 

The  factor  jg"*  in  the  above  inequality  places  such  a  stringent  restriction 
on  step  sizes  that  either  the  increments  in  "z  become,  time  wise,  intol¬ 
erably  small  or  the  increments  in  ^  becomes  too  large  for  acceptable 
truncation  errors. 

2.  The  Explicit  Scheme  of  DuFort  and  Frankel  (D-F  Scheme) 

(A)  Difference  Equation  -  Referring  to  the  sketch  (b)  of  Figure  1,  for 
the  D-F  scheme,  the  unknown  quantity  at  a  grid  point  (  Mt/j  'K/  )  is 

evaluated  in  terms  of  the  known  quantities  at  the  points  (  Ttt!  ), 

(  7L-1  )  and  (  ifL-l ^  K/  ).  Determination  of  the  value  at  the  new 
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point  (  K/)  therefore  involves  data  along  two  preceding  grid  lines, 

i .  e .  the  ^  and  (  )  columns.  The  difference  quotients  for  the 

2  3  4 

D-F  scheme  are  the  first  terms  in  the  following  developments  ’  ’ 

~  ^ -  ^ 

Z  (ZTi)  &  * 

—  .Wfh  Kf(  ~  J  -  - - 

J¥  zc^^)  i 

- 

=  . ' 

The  difference  quotients  used  in  the  D-F  scheme,  which  do  not  involve  the 
point  (>«/,«/)  differ  from  those  of  Eq.  (7)  in  the  E  scheme.  With  Eq.  (17) 
the  differential  Eq.  (1)  can  be  satisfied  approximately  at  the  point  (  tfty ,  rty  ), 
The  following  difference  equation  is  obtained 

where 


.  ,  3/m::z  Czz) 

II 

8/z7^(zx) 

/J,  / - 

II 

"l- 

2  pty 

In  order  to  determine  ^  at  the  unknown  point  (  Tfl-i-l  t  /t'  ),  values  of  ^ 
at  four  points  of  preceding  grid  lines  (  W.,  ft-i-l  ),  (  >»*->  TV  ),  (  ift^  ri’-l  ) 

and  (  Tfi-"/  f  Tty  )  have  to  be  used.  It  may  be  noted  that  if  the  coefficients  in 
the  differential  equation  were  constants,  the  information  at  the  point  (  /«-,  tv) 
would  not  be  required,  since  the  difference  quotients  of  this  scheme  do  not 
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involve  the  point  (  nJ) . 

The  above  difference  Equation  (18)  is  stable  for  all  values  of  (  ) 

and  (  ),  and  the  results  to  be  cited  here  corroborate  this  conclusion. 

Some  further  discussion  on  the  stability  of  the  D-F  scheme  for  the  case 
of  constant  coefficient  has  been  described  by  Forsythe^  and  by  Richtmeyer.  ^ 


(B)  Truncation  Error  and  Convergence 

Assuming  that  IV  is  regular,  the  difference  equation  Eq,  (18)  is 
subject  to  a  truncation  error  of  (AX .  The  term 

0 which  results  from  the  remainder  in  the  second- 

order  difference  quotient  of  Eq.  (17),  reveals  an  undesirable  feature  of  the 
D-F  scheme.  That  is,  the  truncation  error  introduced  by  the  difference 
quotients  cannot  be  reduced  by  considering  arbitrarily  small  AX  and 
alone,  but  rather  it  is  necessary  to  impose  the  condition  AX«A'\p'<C  1  ; 
otherwise,  the  difference  equation  will  not  converge  to  the  parabolic  equation, 
but  to  one  of  a  hyperbolic  type. 

In  order  to  estimate  the  magnitudes  of  the  truncation  error  contributed 
by  various  terms  of  a  difference  equation  similar  to  that  of  Eq.  (11)  may  be 
written  to  include  the  remainders  of  Eq.  (17)  in  the  following  form 


Wz  -  2  f-ivv  -  4  'P  d'w  = 


(19) 


where  and  the  difference  quotients  which  in  this  case  are 

those  defined  by  the  first  terms  of  Eq.  (17).  The  terms  on  the  right-hand 
side  of  Eq.  (19)  constitute  the  truncation  error  of  the  difference  equation 
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Eq.  (18).  As  <? — *-0  ,  and  with  the  aid  of  Eq.  (5)  the  seven  terms  of 
Eq.  (19)  can  be  shown  to  be  in  the  proportion 

ad),  oil)- o(i):0l%]--  o[0^  uo) 

Clearly,  the  presence  of  the  singularity  at  jc  =  0  not  only  magnified  the 
truncation  error  in  the  D-F  scheme  but  also  makes  the  problem  of  conver- 
gence  more  critical. 


3.  .  Implicit  Scheme  (  I  Scheme) 

(A)  Difference  Equation 

In  this  scheme,  the  unknown  quantities  along  a  column  of  grid  points, 
say  the  (  juty-/ )- column,  are  solved  simultaneously  with  known  data  along 
the  fn>  Si'  grid  line  (refer  to  sketch  (c)  of  Fig.  1).  The  difference  quotients 
defining  the  I- scheme  are  constructed  according  to  the  first  terms  in  the 

2  9 

following  expansion  * 


^  - L/^^)  U/_,  -f . . 

^  - 

^'ip'  Z(A<P)  ^ 

^IV  =:  WmH.M  ~Z.  lA/^  -f - 


>-  (21) 


With,  the  difference  quotients  of  Eq.  (21),  the  differential  equation,  Eq.  (1), 
is  reduced  to  a  system  of  linear  algebraic  equations  with  the  in  the 
(  )-column. 


For  numerical  consideration,  one  must  multiply  the  fifth  and  the  sixth  terms 
in  Eq.  (20)  by  a  factor  of  1  /30.  Note  also  that  >  A‘ip‘ 
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(22) 


^Ot)  Whh,k  ^Cft) 


where 


4</zi;^(^x) 


/ 


S^2W,^.  A.  Cax) 


4)l2Wm,l/'Ax) 


The  integer  fV  is  taken  from  1  to  N  for  each  nv  ,  with  XV  =  I  corresponding 
to  the  first  grid  point  from  the  "cone  surface"  and  n»A/  corresponding  to  the 
boundary  point  at  the  shock.  Or,  in  an  expansion  form 


j-  Z  ’f’ j(^)W)iltly3 


1^,/  1 

=  Wh,Z 


where  W^fijO  denote  values  to  be  determined  by  the  inner 

and  outer  boundary  conditions  respectively.  These  equations,  together  with 
the  two  boundary  conditions,  form  a  system  which  suffices  for  the  determin¬ 
ation  of  the  (  /\/+ 1  )  unknowns.  The  matrix  of  this  system  is  of  the  tridiagonal 
type,  and  the  solution  can  thus  be  obtained  by  following  standard  procedures. 
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2 

e.g.  Gauss'  elimination  method. 

(B)  Truncation  Error 

The  choice  of  the  mesh  spacings  for  AS  and  in  this  scheme  is 

2  3 

not  subject  to  the  stability  and  the  convergence  requirements.  ’  Generally, 
in  the  I-scheme,  the  difference  equation,  Eq.  (22),  is  subject  to  the  trunca¬ 
tion  errors  of  the  order  •  The  error  of  arises  from 

the  use  of  the  backward-difference  quotient  for  /dS  as  well  as  from  re¬ 
placing  y  Z  Wifui  li  y  2 n.  •  Note  that  the  differential  equation  is 

satisfied  in  this  scheme  at  the  station  (  !  ).  To  examine  more  thoroughly 

the  manner  by  which  the  singularity  at  ^  =  0  affects  the  truncation  error, 
one  may  again  express  the  difference  equations,  Eq.  (22),  in  a  form  similar 
to  that  of  Eq.  (11)  and  (19). 

where  the  terms  on  the  right-hand  side  again  represent  the  remainder  of  the 
difference  equation.  The  same  estimate  of  each  term  in  Eq.  (11)  holds  for 
the  corresponding  term  in  Eq.  (23)  except  the  last  (additional)  term,  which 
arises  from  replacing  y  Z Wm*i, ^  order  of 

.  As  tC - >■  O  *  th®  seven  terms  of  Eq.  (23)  are  thus  in  the 

proportion 

oa)  *  oa)  --od)' o[f§-]  ^ .*  0  ’  o  (24) 

Clearly  the  accuracy  of  this  scheme  is  still  greatly  affected  by  the  singular 
point.  The  problem,  however,  is  far  less  critical  than  encountered  in  the  E 


(23) 
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and  D-F  schemes. 

4.  Boundary  Conditions 

The  conditions  at  the  "cone  surface"  and  at  the  shock  are  the  inner  and 
outer  boundary  conditions,  respectively.  The  difference  equation  for  the 
inner  boundary  condition  is  simply 

(25) 

In  order  to  determine  IV  at  the  outer  boundary  with  an  accuracy  comparable 
to  that  obtained  with  the  difference  equations,  the  derivative  a^t  a  boun¬ 

dary  point  will  be  expressed  in  terms  of  iV  at  the  neighboring  points  (refer 
to  sketch  (d)  in  Fig.  1)  as 

=  .  >3  Wg^U  ^  ’  ^  (26) 

N 

Notice  that,  because  of  the  term  y2(V '  >  the  boundary  condition  Eq.  (3)  is 
nonlinear.  Hence,  one  would  have  to  solve  a  quadratic  algebraic  equation 
for  the  boundary  value  of  pV  •  To  simplify  the  numerical  computation,  one 
may  replace  V ^  A/  ^  ^  N  •  which  introduces  only  an  error 

of  and  is  consistent  with  the  accuracy  of  all  three  schemes.  The 

outer  boundary  conditions  can  therefore  be  satisfied 

In  the  implicit  scheme,  Eq.  (27)  is  used  to  eliminate  N 

equation  of  the  system  of  Eq.  (22). 


Again,  in  order  to  compare  the  relative  magnitudes  of  the  terms  in  Eq. 
(24)  in  a  strict  numerical  sense,  one  must  multiply  both  ^  and 

^  factor  of  1/30. 


( 
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V.  DISCUSSION  OF  RESULTS 


In  the  preceding  sections,  three  difference  schemes  have  been  des¬ 
cribed  for  the  cone  problem.  Excluding  the  region  near  the  singularity 
X  =  0,  the  truncation  error  in  the  difference  approximation  belongs  to 
0  [(^X) (Alp/]  for  the  E  and  I  schemes,  and  to  0 for  the 
D-F  scheme.  By  making  AX  «  A'p'«l  ,  not  only  can  the  truncation  errors 
in  all  the  schemes  be  reduced,  but  also  the  stability  and  convergence  require¬ 
ments  for  the  E  and  D-F  schemes  can  be  satisfied. 

The  singularity  at  ^  =  0  makes  the  accuracy  of  the  solution,  or 
alternatively  the  efficiency  of  all  three  schemes,  very  poor,  as  the  cone  apex 
is  approached.  The  singularity  also  makes  the  stability  and  convergence  pro¬ 
blems  of  the  E  and  D-F  schemes  more  critical.  This  poses  a  serious  problem 
in  practice,  because  application  of  any  of  these  schemes  over  a  range  of  ^  , 
which  begins  at  a  large  ^  to  avoid  the  adverse  effect  of  the  singularity, 
cannot  avoid  encountering  a  difficulty  of  another  kind.  For,  in  order  to  provide 
accurate  initial  data  at  a  station  far  removed  from  "X  =  0,  one  has  to  use 
a  large  number  of  high-order  terms  in  the  development  of  Eq.  (5).  In  Eq,  (5), 
a  four-term  expansion  has  been  obtained  for  the  purpose  of  the  present  study. 
But  in  other  more  general  problems,  the  task  of  determining  the  coefficients 
of  the  higher-order  terms  would  be  too  burdensome. 

In  the  following  comparison  of  the  three  schemes,  examination  will  be 
made  first  over  the  range  O.Z^  X  J.  ,  and  then  over  a  more  critical  range 
,0!  X  ^  0.2.  The  skin-friction  and  heat-transfer  coefficients,  as  well 

as  the  velocity  profiles  obtained  by  the  I- scheme  will  be  presented  over  a 
wide  range:  0.  01  <  ^ 
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5. 1  Niutierical  Integration  Beginning  at  =0.  20 
From  the  Crocco  relation  Eq.  (4), 

Cm  =  Cp/ Z  COS^ 


(28) 


where  Cm  and  Cp  are  the  surface  heat-transfer  and  skin-friction  coefficients, 
defined  as 


Cm  : 

Cf 


Ua,  z  a  ip 


(29) 


Typical  results  obtained  from  the  three  methods  are  presented  in  Figs.  2  and 
3.  In  Fig,  2,  Cp  and  Cp  obtained  by  the  E,  I  and  D-F  schemes  over  the 
range  0,  2  ^  ^  ^  1.0  are  compared.  In  order  to  provide  the  initial  data  at 
^=0.2  with  an  error  no  more  than  1%,  all  terms  in  the  four-term  expan¬ 
sion  in  5?  given  in  Eq.  (5)  have  to  be  used.  The  calculations  were  performed 
by  an  IBM  704  digital  computer,  carrying  eight  significant  figures. 

It  was  found  that  the  choice  of  -  1  /lO  is  sufficient  for  all  three 

schemes.  However,  the  step-size  had  to  be  varied  with  different  schemes 

to  accomodate  the  stability  and  convergence  requirements  as  well  as  the  singu- 

-5 

larity  affect  on  truncation  errors.  In  the  E  scheme,  the  increment  5x1 0 

is  required  in  order  to  fulfill  the  stability  condition  at  ^  =  .20.  A  smaller 

-5 

step-size  Ax  =  2.5  x  10  was  chosen  in  the  calculation  because  instability 

—  -5 

has  actually  occurred  with  Ax  =  5x10  .  Two  runs  were  made  for  each  of 

the  I  and  D-F  schemes,  first  with  Ax  -  1/100,  =1/10  and  the  second 

with  Ax  =  1/1000,  A^-  1/10.  The  good  agreement  among  results  of  the 

second  runs  of  I  and  D-F  schemes  as  well  as  the  E  scheme  provides  a  prelim¬ 
inary  indication  of  the  adequate  accuracy  of  these  schemes.  The  accuracies 
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I 


are  less  satisfactory  for  the  first  runs  of  the  I  and  D-F  schemes  with 

=  1/100,  particularly  the  latter.  The  total  computation  times  (includ¬ 
ing  print-out  at  the  .05-intervals  of  x  )  are  approximately  45  minutes 
for  the  E  scheme,  20  seconds  for  the  first  runs  of  the  I  and  D-F  schemes, 
and  100  seconds  for  the  second  runs  of  the  I  and  D-F  schemes. 

Notice  that  the  running  time  for  both  I  and  D-F  schemes  are  about 
the  same;  nevertheless,  with  the  same  set  of  increments,  the  D-F  scheme 
is  less  accurate  as  shown  in  Fig.  2. 

5.  2  Numerical  Integration  Beginning  at  ^  =  .  01 

As  a  more  crucial  test  of  the  usefulness  of  the  methods,  calculations 
were  carried  out  and  compared  over  the  range  of  a  smaller  x  (0- 1  ^  ^  ^  0.  3). 
The  initial  data  can  be  adequately  provided  by  the  first  two-terms  in  the 
development  of  £q.  (5)  but  the  problems  of  the  singularity  tends  to  be  more 
critical.  In  fact,  the  computation  by  IBM  704  for  the  E  scheme  becomes 
virtually  impossible  under  the  stringent  stability  requirement  over  this  range. 
Thus,  only  the  results  from  I  and  D-F  schemes  are  considered.  The  two 
schemes  are  first  carried  out  with  the  increments  AX  =  1/100,  and  1/10, 
then  with  =  1/1000  and  =  1/10.  As  shown  in  Fig.  3,  the  departure 
of  the  solution  by  the  D-F  scheme  is  greatly  amplified  at  small  7  ,  where 
the  singularity  makes  satisfaction  of  the  convergence  requirement  more 
difficult.  Even  in  the  second  run  with  a  more  refined  step- size  in  ^  ,  a 
small  but  appreciable  error  of  about  4%  still  exists  for  the  D-F  scheme. 

The  I- scheme  in  the  second  run  with  1/1000  and  =  1/10  yields 

good  accuracy,  agreeing  with  the  four-term  analytic  development  of  Eq.  (5) 
throughout  the  entire  range  of  .  01  <  .  3.  The  computation  time 
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(including  print-out)  of  the  runs  with  -  1/1000  and  aW  -  1/10  was 
about  110  seconds  for  each  scheme. 

5.  3  Discussion  on  the  Efficiency  of  the  Three  Schemes 

In  summarizing  the  results  presented  in  Figs.  2  and  3,  the  £  scheme 
is  definitely  less  efficient  and  not  suitable  for  the  treatment  of  singular  dif¬ 
ferential  equations  such  as  Eq.  (1).  With  the  same  intervals  of  AX  and  A^ , 
the  computation  time  for  the  E  scheme  is  45  minutes,  but  only  20  seconds  for 
the  I  and  D-F  schemes.  The  examples  show  that  the  singularity  at  X  -  Q 
affects  critically  the  accuracy  and  efficiency  of  both  the  I  and  the  D-F  schemes, 
especially  if  higher-order  terms  in  the  development  about  the  singular  point 
are  not  available.  The  I  scheme,  which  is  subject  to  neither  the  stability  nor 
the  convergence  requirements,  proves  in  this  instance  to  be  somewhat  more 
reliable  than  the  D-F  scheme.  This  relative  merit  of  the  I  scheme  is  com¬ 
prehensible  in  the  light  of  the  discussion  in  Section  IV.  One  may  recall  that 
the  largest  truncation  error  in  the  I  scheme  is  of  0(~^)  whereas  the  largest 
in  the  D-F  scheme  is  of  - 

It  is  of  interest  to  note  that  the  truncation  errors  associated  with  A'^ 
were  found  to  be  very  small.  Calculations  have  been  repeated  with  A^ 
changed  from  1/10  to  1/50  without  appreciable  difference  in  the  solution. 

The  comparative  insensitivity  of  the  solutions  to  the  interval  in  both  the 

I  and  D-F  schemes  is  also  understandable  from  the  previous  discussion, 

2 

that  there  is  always  a  small  numerical  factor  associated  with  the  (  A^  ) 
terms  in  the  remainders  of  the  difference  equation,  which  compensates 
partly  for  the  adverse  effects  of  the  singularity.  This  insensitivity  is,  of 
course,  a  property  peculiar  to  the  specific  problem  considered.  However, 
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in  a  more  general  case,  especially  if  the  pressure  gradients  do  not  vanish, 

one  may  have  to  refine  not  only  the  increment  in  %  but  also  the  step- size 
_  _  >)< 

of  as  X - ^  Q  .  In  that  case,  the  freedom  from  the  convergence 

requirement  makes  the  I  scheme  even  more  superior. 

5,4  Discussion  of  the  Numerical  Solutions 

In  order  to  demonstrate  the  accuracy  of  the  method  the  numerical 
solution  to  Eqs.  (1)  -  (3)  is  carried  out  for  hypersonic  flow  over  a  cone 
under  the  I  scheme  from  ^  =  0.  01  to  %  =  10.  The  computation  was  made 
in  two  stages  of  different  step  sizes  in  ^  .  In  the  first  stage,  the  integra¬ 
tion  by  the  I  scheme,  with  AX  -  1/1000,  A^=  1/10,  was  carried  out  from 

"it  =  0.01  to  7  =  0.2  with  the  initial  data  provided  by  the  four-term  expan¬ 
sion  of  Eq.  (5).  The  integration  by  I  scheme  was  continued  in  the  second 
stage  from  Jc  -  0.2  through  10  with  a  larger  step- size  in  ^  ,  (  AX  =  1/100, 
A^  =  1  /lO).  The  total  computation  time  {including  print-out  at  0.  05  interval 
of  7  )  was  about  two  minutes,  with  38  seconds  for  the  first  stage  and  90 
seconds  for  the  second  stage.  The  results  for  skin  friction,  surface  heat- 
transfer  rate,  and  velocity  profiles  are  presented  in  Figs.  4,  5,  and  6. 

The  ratios  of  and  given  in  Fig.  4  are  in  good  agree¬ 

ment  with  the  analytic  development  obtained  for  the  region  of  X<L  ,  and 
approach  the  correct  values  of  the  classical  boundary -layer  theory  for 

5r - ►■00  .  The  velocity  profiles  of  the  shock  layer  at  various  distances 

from  the  apex  are  given  in  Fig.  5.  At  low  value  of  jt  .  the  distribution  is 
_________ 

In  such  a  general  case,  the  correction  for  the  singularity  associated  with 
^  =  0  will  be  no  longer  as  small  as  represented  by  Eq.  (10),  and  cannot 
be  neglected.  However,  as  pointed  out  previously,  this  singvilarity  can  be 
eliminated  by  using  the  variable  VC  and  instead  of  W  and  ^  .  The  new 
variables  will  introduce  a  factor  ofl/f^to  the  remainder  of  the  difference 
equations  and  make  the  use  of  the  I  scheme  (rather  than  other  schemes)  more 
desirable. 
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% 


linear  in  corresponding  to  a  uniform  shear  flow.  As  X  increases, 

the  velocity  gradient  increases  near  the  surface  and  decreases  near  the 
shock,  with  the  velocity  in  the  outer  portion  approaching  to  the  tangential 
component  of  the  free- stream  velocity. 

It  is  rather  interesting  to  observe  that,  even  at  7  as  low  as  3, 
the  flow  field  can  already  be  represented  very  closely  by  the  classical 
boundary-layer  limit,  in  spite  of  the  fact  that  the  velocity  at  the  outer  edge 
is  still  quite  far  from  the  inviscid  value  and  that  the  viscous  layer  is  still 
a  major  fraction  of  the  shock  layer.  To  clearly  demonstrate  this  observation, 
the  velocity  profiles  at  various  X  stations  are  correlated  in  terms  of  the 
similarity  variable 

^  _  rjEZ 

V  “  I'tr i  {ZW 

and  presented  in  Fig.  6.  The  good  correlation  and  excellent  agreement  with 
the  Blasius  profile  provide  another  check  on  the  validity  and  accuracy  of  the 
method. 

It  should  be  noted  that,  in  order  to  apply  the  solution  obtained,  one 

must  relate  the  variable  %  with  the  Reynolds  number  (or  its  reciprocal, 

the  Knudsen  number)  based  on  the  free- stream  flow  condition.  This  involves 

a  choice  of  the  reference  temperature  7#.  because  of  the  linear  viscosity- 

1  15 

temperature  relation  used,  ’  which  in  this  case  may  be  simply  taken  as 


I 

{  * 
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VI.  CONCLUSION 


The  relative  merit  of  three  difference  schemes  for  integrating  the 
viscous  shock-layer  equation  has  been  investigated  with  specific  reference 
to  the  problem  of  low  density  hypersonic  flow  over  a  cone.  The  difference 
quotients  and  equations  defining  each  of  the  schemes  have  been  derived  and 
the  nature  of  the  singularities  on  the  solution  and  their  effects  on  the  stability 
and  convergence  characteristics  of  these  schemes  were  studied. 

There  are  two  singularities  in  the  solution.  One  is  associated  with 
the  leading-edge  (  =  0)  and  the  other  with  the  cone  surface  {  '\p'=  0), 

Close  examination  shows  that  the  latter  has  little  effect  on  the  results.  The 
leading-edge  singularity  on  the  other  hand  critically  affects  the  accuracy  and 
efficiency  of  all  schemes,  particularly  the  stability  characteristics  of  the 
explicit  scheme.  Because  of  this  effect,  one  must  either  increase  the  nvimber 
of  the  higher-order  terms  in  the  power- series  expansion  around  the  leading 
edge  (so  as  to  maintain  the  accuracy  of  the  initial  data  at  a  comparatively 
large  value  of  ^  ),  or  further  refine  the  step- size  in  ^  (so  as  to  maintain 
the  accuracy  of  the  difference  equation  at  a  comparatively  small  value  of  /T  ). 

The  three  schemes  applied  to  the  cone  problem  were  carried  out  with 
the  aid  of  the  IBM  704  digital  computer.  The  results  indicate  that,  because  of 
the  leading-edge  singularity,  the  stability  requirement  becomes  too  stringent 
for  the  explicit  scheme  to  be  useful.  While  the  difference  in  efficiency  of  the 
implicit  and  OuFort-Frankel  schemes  is  actually  not  large  for  this  example, 
the  implicit  scheme,  which  is  subject  to  neither  the  stability  nor  the  conver¬ 
gence  requirements,  is  seen  to  be  more  reliable.  It  is  noted,  furthermore, 
that  in  the  more  general  problems  with  non-vanishing  pressure  gradients,  the 
implicit  scheme  should  be  more  advantageous. 
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(a)  THE  EXPLICIT  SCHBIE  (b)  THE  EXPLICIT  SCHEME  OF 

DUFORT-FRANKEL  TYPE 


X  UMKHONM  POIMT 
•  KHONH  POIMT 


(c)  THE  IMPLICIT  SCHEME  (d)  OIFFEREMCE  ^TIMEMT  FOR 

THE  BOUHDARV  POIMT 


Figure  I  SKETCH  FOR  COHSTRUCTINQ  DIFFERENCE  QUOTIENTS 
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Figure  3  SURFACE  HEAT-TRANSFER  AND  SKIN-FRICTION  COEFFICIENTS  AS  DETERMINED 
BY  THE  IMPLICIT  AND  OUFORT-FRANXa  SCHEMES 
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Figure  4  SURFACE  HEAT  TRANSFER  AND  SKIN  FRICTION  AT  VARIOUS  DEGREES  OF  GAS  RAREFACTION 


